{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Comparison of the performance of the Circuit-constraint\n",
    "\n",
    "There are multiple ways to model the Traveling Salesman Problem (TSP) in CP-SAT.\n",
    "In this notebook, we compare the performance of some."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "TIMELIMIT_S = 10\n",
    "STEP_SIZE = 10"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Generate random instances for a benchmark\n",
    "\n",
    "To compare the performance of the different models, we need to generate random instances.\n",
    "Because purely random graphs are not very interesting, we use points in the plane and the distance between them as edge weights.\n",
    "This gives the instances a more realistic structure.\n",
    "Without structure, solvers like CP-SAT and Gurobi cannot be smart and the comparison would be boring."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "import networkx as nx\n",
    "import random\n",
    "import itertools\n",
    "import typing\n",
    "\n",
    "\n",
    "# Generate a random graph with a specified number of nodes.\n",
    "def generate_random_euclidean_graph(\n",
    "    n: int,\n",
    ") -> typing.Tuple[nx.Graph, typing.List[typing.Tuple[int, int]]]:\n",
    "    \"\"\"\n",
    "    Generate a random graph with a specified number of nodes.\n",
    "    The nodes will be integers 0, 1, ..., n-1.\n",
    "    The graph will be a complete graph with Euclidean distance as edge weight.\n",
    "    \"\"\"\n",
    "    points = [(random.randint(0, 10_000), random.randint(0, 10_000)) for _ in range(n)]\n",
    "    G = nx.Graph()\n",
    "    G.add_nodes_from(range(n))\n",
    "    # Complete graph with Euclidean distance as edge weight.\n",
    "    for i, j in itertools.combinations(range(n), 2):\n",
    "        p = points[i]\n",
    "        q = points[j]\n",
    "        dist = round((p[0] - q[0]) ** 2 + (p[1] - q[1]) ** 2)\n",
    "        G.add_edge(i, j, weight=dist)\n",
    "    return G, points"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Models"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Model 1: AddCircuit-constraint"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "from ortools.sat.python import cp_model\n",
    "\n",
    "\n",
    "def solve_with_add_circuit(G: nx.Graph) -> typing.Tuple[bool, float]:\n",
    "    \"\"\"Solve TSP with add_circuit.\"\"\"\n",
    "    model = cp_model.CpModel()\n",
    "\n",
    "    # Variables\n",
    "    edge_vars = dict()\n",
    "    for u, v in G.edges:\n",
    "        edge_vars[u, v] = model.new_bool_var(f\"edge_{u}_{v}\")\n",
    "        edge_vars[v, u] = model.new_bool_var(f\"edge_{v}_{u}\")\n",
    "\n",
    "    # Constraints\n",
    "    # Because the nodes in the graph a indices 0, 1, ..., n-1, we can use the\n",
    "    # indices directly in the constraints. Otherwise, we would have to use\n",
    "    # a mapping from the nodes to indices.\n",
    "    circuit = [(u, v, x) for (u, v), x in edge_vars.items()]\n",
    "    model.add_circuit(circuit)\n",
    "\n",
    "    # Objective\n",
    "    model.minimize(sum(x * G[u][v][\"weight\"] for (u, v), x in edge_vars.items()))\n",
    "\n",
    "    solver = cp_model.CpSolver()\n",
    "    solver.parameters.max_time_in_seconds = TIMELIMIT_S\n",
    "    # solver.parameters.log_search_progress = True\n",
    "    status = solver.solve(model)\n",
    "    return status == cp_model.OPTIMAL, solver.objective_value\n",
    "\n",
    "\n",
    "assert solve_with_add_circuit(generate_random_euclidean_graph(10)[0])[0], (\n",
    "    \"Small graph should be solvable.\"\n",
    ")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Model 2: Miller Tucker Zemlin"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "def solve_with_miller_tucker_zemlin(G: nx.Graph) -> typing.Tuple[bool, float]:\n",
    "    model = cp_model.CpModel()\n",
    "\n",
    "    # Variables\n",
    "    edge_vars = dict()\n",
    "    for u, v in G.edges:\n",
    "        edge_vars[u, v] = model.new_bool_var(f\"edge_{u}_{v}\")\n",
    "        edge_vars[v, u] = model.new_bool_var(f\"edge_{v}_{u}\")\n",
    "\n",
    "    depth_vars = dict()\n",
    "    for u in G.nodes:\n",
    "        depth_vars[u] = model.new_int_var(0, G.number_of_nodes() - 1, f\"depth_{u}\")\n",
    "\n",
    "    # Constraints\n",
    "    # Every node has exactly one incoming and one outgoing edge.\n",
    "    for v in G.nodes:\n",
    "        model.add_exactly_one([edge_vars[u, v] for u in G.neighbors(v)])\n",
    "        model.add_exactly_one([edge_vars[v, u] for u in G.neighbors(v)])\n",
    "\n",
    "    # Use depth variables to prohibit subtours.\n",
    "    model.add(depth_vars[0] == 0)  # fix the root node to depth 0\n",
    "    for v, w in G.edges:\n",
    "        if w != 0:  # The root node is special.\n",
    "            model.add(depth_vars[v] + 1 == depth_vars[w]).only_enforce_if(\n",
    "                edge_vars[v, w]\n",
    "            )\n",
    "        if v != 0:\n",
    "            model.add(depth_vars[w] + 1 == depth_vars[v]).only_enforce_if(\n",
    "                edge_vars[w, v]\n",
    "            )\n",
    "\n",
    "    # Objective\n",
    "    model.minimize(sum(x * G[u][v][\"weight\"] for (u, v), x in edge_vars.items()))\n",
    "\n",
    "    solver = cp_model.CpSolver()\n",
    "    solver.parameters.max_time_in_seconds = TIMELIMIT_S\n",
    "    # solver.parameters.log_search_progress = True\n",
    "    status = solver.solve(model)\n",
    "    return status == cp_model.OPTIMAL, solver.objective_value\n",
    "\n",
    "\n",
    "assert solve_with_miller_tucker_zemlin(generate_random_euclidean_graph(10)[0])[0], (\n",
    "    \"Small graph should be solvable.\"\n",
    ")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "n=10\n",
      "  add_circuit solved it with objective 185396654.0.\n",
      "  miller_tucker_zemlin solved it with objective 185396654.0.\n",
      "n=20\n",
      "  add_circuit solved it with objective 109773872.0.\n",
      "  miller_tucker_zemlin solved it with objective 109773872.0.\n",
      "n=30\n",
      "  add_circuit solved it with objective 83156520.0.\n",
      "  miller_tucker_zemlin solved it with objective 83156520.0.\n",
      "n=40\n",
      "  add_circuit solved it with objective 85448930.0.\n",
      "  miller_tucker_zemlin solved it with objective 85448930.0.\n",
      "n=50\n",
      "  add_circuit solved it with objective 106851566.0.\n",
      "  miller_tucker_zemlin dropped out.\n",
      "n=60\n",
      "  add_circuit solved it with objective 74616450.0.\n",
      "n=70\n",
      "  add_circuit solved it with objective 88725432.0.\n",
      "n=80\n",
      "  add_circuit solved it with objective 77799112.0.\n",
      "n=90\n",
      "  add_circuit dropped out.\n"
     ]
    }
   ],
   "source": [
    "challengers = {\n",
    "    \"add_circuit\": solve_with_add_circuit,\n",
    "    \"miller_tucker_zemlin\": solve_with_miller_tucker_zemlin,\n",
    "}\n",
    "n = STEP_SIZE\n",
    "while challengers:\n",
    "    print(f\"n={n}\")\n",
    "    G, points = generate_random_euclidean_graph(n)\n",
    "    dropouts = []\n",
    "    for name, solver in challengers.items():\n",
    "        solved, obj = solver(G)\n",
    "        if not solved:\n",
    "            dropouts.append(name)\n",
    "        else:\n",
    "            print(f\"  {name} solved it with objective {obj}.\")\n",
    "    for name in dropouts:\n",
    "        print(f\"  {name} dropped out.\")\n",
    "        del challengers[name]\n",
    "    n += STEP_SIZE"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Model 3: Iterative Dantzig Model\n",
    "\n",
    "Usually, you would use callbacks to add constraints dynamically.\n",
    "However, this is not possible in the CP-SAT solver.\n",
    "Because adding all constraints from the beginning is prohibitive, we add the constraints after each round, and solve the model again."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "import time\n",
    "\n",
    "\n",
    "class SubtourCallback(cp_model.CpSolverSolutionCallback):\n",
    "    \"\"\"\n",
    "    Callback to detect subtours during the search.\n",
    "    \"\"\"\n",
    "\n",
    "    def __init__(self, edge_vars: typing.Dict, early_abort=False):\n",
    "        \"\"\"\n",
    "        edge_vars: a dictionary mapping edges to variables.\n",
    "        early_abort: if True, the search will be aborted as soon as a subtour is detected. This is closer to lazy constraints, but may waste the previous work.\n",
    "        \"\"\"\n",
    "        super().__init__()\n",
    "        self.edge_vars = edge_vars\n",
    "        self.subtours = []\n",
    "        self.early_abort = early_abort\n",
    "\n",
    "    def on_solution_callback(self):\n",
    "        connected_components = list(\n",
    "            nx.connected_components(\n",
    "                nx.Graph(\n",
    "                    (u, v) for (u, v), x in self.edge_vars.items() if self.value(x)\n",
    "                )\n",
    "            )\n",
    "        )\n",
    "        if len(connected_components) > 1:\n",
    "            self.subtours += connected_components\n",
    "            if self.early_abort:\n",
    "                self.stop_search()\n",
    "\n",
    "    def has_subtours(self):\n",
    "        return len(self.subtours) > 0\n",
    "\n",
    "    def reset(self):\n",
    "        self.subtours = []\n",
    "\n",
    "\n",
    "def solve_with_dantzig(G: nx.Graph) -> typing.Tuple[bool, float]:\n",
    "    model = cp_model.CpModel()\n",
    "\n",
    "    # Variables\n",
    "    edge_vars = dict()\n",
    "    for u, v in G.edges:\n",
    "        # symmetric variables\n",
    "        x = model.new_bool_var(f\"edge_{u}_{v}\")\n",
    "        edge_vars[u, v] = x\n",
    "        edge_vars[v, u] = x\n",
    "\n",
    "    # Constraints\n",
    "    # Every nodes has two incident edges.\n",
    "    for v in G.nodes:\n",
    "        model.add(sum(edge_vars[u, v] for u in G.neighbors(v)) == 2)\n",
    "\n",
    "    # Objective\n",
    "    tour_cost = sum(G.edges[u, v][\"weight\"] * edge_vars[(u, v)] for u, v in G.edges)\n",
    "    model.minimize(tour_cost)\n",
    "\n",
    "    solver = cp_model.CpSolver()\n",
    "\n",
    "    start_time = time.time()\n",
    "\n",
    "    def remaining_time():\n",
    "        return TIMELIMIT_S - (time.time() - start_time)\n",
    "\n",
    "    solver.parameters.max_time_in_seconds = remaining_time()\n",
    "    callback = SubtourCallback(edge_vars, early_abort=False)\n",
    "    # solver.parameters.log_search_progress = True\n",
    "    status = solver.solve(model, callback)\n",
    "\n",
    "    # The following part is more complex. Here we repeatedly add constraints\n",
    "    # and solve the model again, until we find a solution without subtours.\n",
    "    # TODO: Use the previous solutions as hints for the next solve.\n",
    "    while callback.has_subtours():\n",
    "        subtours = callback.subtours\n",
    "        if len(subtours) > 1:\n",
    "            for comp in subtours:\n",
    "                outgoing_edges = sum(\n",
    "                    edge_vars[u, v]\n",
    "                    for u in comp\n",
    "                    for v in G.neighbors(u)\n",
    "                    if v not in comp\n",
    "                )\n",
    "                model.add(outgoing_edges >= 2)\n",
    "            callback.reset()\n",
    "            tour_cost = sum(\n",
    "                G.edges[u, v][\"weight\"] * edge_vars[(u, v)] for u, v in G.edges\n",
    "            )\n",
    "            model.add(\n",
    "                tour_cost >= int(solver.best_objective_bound)\n",
    "            )  # help with lower bound\n",
    "            solver.parameters.max_time_in_seconds = remaining_time()\n",
    "            if remaining_time() <= 0:\n",
    "                return False, solver.objective_value\n",
    "            status = solver.solve(model, callback)\n",
    "        else:\n",
    "            break\n",
    "    return status == cp_model.OPTIMAL, solver.objective_value\n",
    "\n",
    "\n",
    "assert solve_with_dantzig(generate_random_euclidean_graph(20)[0])[0], (\n",
    "    \"Small graph should be solvable.\"\n",
    ")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "n=10\n",
      "  add_circuit solved it with objective 97715906.0.\n",
      "  miller_tucker_zemlin solved it with objective 97715906.0.\n",
      "  dantzig solved it with objective 97715906.0.\n",
      "n=20\n",
      "  add_circuit solved it with objective 86770784.0.\n",
      "  miller_tucker_zemlin solved it with objective 86770784.0.\n",
      "  dantzig solved it with objective 86770784.0.\n",
      "n=30\n",
      "  add_circuit solved it with objective 103530600.0.\n",
      "  miller_tucker_zemlin solved it with objective 103530600.0.\n",
      "  dantzig solved it with objective 103530600.0.\n",
      "n=40\n",
      "  add_circuit solved it with objective 76480910.0.\n",
      "  miller_tucker_zemlin solved it with objective 76480910.0.\n",
      "  dantzig solved it with objective 76480910.0.\n",
      "n=50\n",
      "  add_circuit solved it with objective 77015310.0.\n",
      "  miller_tucker_zemlin solved it with objective 77015310.0.\n",
      "  dantzig solved it with objective 77015310.0.\n",
      "n=60\n",
      "  add_circuit solved it with objective 88887270.0.\n",
      "  dantzig solved it with objective 88887270.0.\n",
      "  miller_tucker_zemlin dropped out.\n",
      "n=70\n",
      "  add_circuit solved it with objective 82404140.0.\n",
      "  dantzig dropped out.\n",
      "n=80\n",
      "  add_circuit solved it with objective 76343838.0.\n",
      "n=90\n",
      "  add_circuit solved it with objective 74842530.0.\n",
      "n=100\n",
      "  add_circuit solved it with objective 85180224.0.\n",
      "n=110\n",
      "  add_circuit solved it with objective 74507340.0.\n",
      "n=120\n",
      "  add_circuit dropped out.\n"
     ]
    }
   ],
   "source": [
    "challengers = {\n",
    "    \"add_circuit\": solve_with_add_circuit,\n",
    "    \"miller_tucker_zemlin\": solve_with_miller_tucker_zemlin,\n",
    "    \"dantzig\": solve_with_dantzig,\n",
    "}\n",
    "n = STEP_SIZE\n",
    "while challengers:\n",
    "    print(f\"n={n}\")\n",
    "    G, points = generate_random_euclidean_graph(n)\n",
    "    dropouts = []\n",
    "    for name, solver in challengers.items():\n",
    "        solved, obj = solver(G)\n",
    "        if not solved:\n",
    "            dropouts.append(name)\n",
    "        else:\n",
    "            print(f\"  {name} solved it with objective {obj}.\")\n",
    "    for name in dropouts:\n",
    "        print(f\"  {name} dropped out.\")\n",
    "        del challengers[name]\n",
    "    n += STEP_SIZE"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Model 4: Using a classical MIP-Solver (Gurobi)\n",
    "\n",
    "This is the classical approach to solve the TSP, as the linear programming relaxation is already a good approximation to the solution."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [],
   "source": [
    "import gurobipy as gp\n",
    "\n",
    "\n",
    "def solve_with_gurobi(G: nx.Graph) -> typing.Tuple[bool, float]:\n",
    "    model = gp.Model()\n",
    "    model.Params.LogToConsole = 0\n",
    "    model.Params.TimeLimit = TIMELIMIT_S\n",
    "    model.Params.lazyConstraints = 1\n",
    "\n",
    "    # Variables\n",
    "    edge_vars = dict()\n",
    "    for u, v in G.edges:\n",
    "        # symmetric variables\n",
    "        x = model.addVar(vtype=gp.GRB.BINARY, name=f\"edge_{u}_{v}\")\n",
    "        edge_vars[u, v] = x\n",
    "        edge_vars[v, u] = x\n",
    "\n",
    "    def gurobi_subtour_callback(model, where):\n",
    "        if where == gp.GRB.Callback.MIPSOL:\n",
    "            connected_components = list(\n",
    "                nx.connected_components(\n",
    "                    nx.Graph(\n",
    "                        (u, v)\n",
    "                        for (u, v), x in edge_vars.items()\n",
    "                        if model.cbGetSolution(x) > 0.5\n",
    "                    )\n",
    "                )\n",
    "            )\n",
    "            if len(connected_components) > 1:\n",
    "                for comp in connected_components:\n",
    "                    outgoing_edges = sum(\n",
    "                        edge_vars[u, v]\n",
    "                        for u in comp\n",
    "                        for v in G.neighbors(u)\n",
    "                        if v not in comp\n",
    "                    )\n",
    "                    model.cbLazy(outgoing_edges >= 2)\n",
    "\n",
    "    # Constraints\n",
    "    # Every nodes has two incident edges.\n",
    "    for v in G.nodes:\n",
    "        model.addConstr(sum(edge_vars[u, v] for u in G.neighbors(v)) == 2)\n",
    "\n",
    "    # Objective\n",
    "    tour_cost = sum(G.edges[u, v][\"weight\"] * edge_vars[(u, v)] for u, v in G.edges)\n",
    "    model.setObjective(tour_cost, gp.GRB.MINIMIZE)\n",
    "\n",
    "    model.optimize(gurobi_subtour_callback)\n",
    "    return model.status == gp.GRB.OPTIMAL, model.objVal\n",
    "\n",
    "\n",
    "assert solve_with_gurobi(generate_random_euclidean_graph(20)[0]), (\n",
    "    \"Small graph should be solvable.\"\n",
    ")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "n=10\n",
      "  add_circuit solved it with objective 100129810.0.\n",
      "  miller_tucker_zemlin solved it with objective 100129810.0.\n",
      "  dantzig solved it with objective 100129810.0.\n",
      "  gurobi solved it with objective 100129810.0.\n",
      "n=20\n",
      "  add_circuit solved it with objective 125047848.0.\n",
      "  miller_tucker_zemlin solved it with objective 125047848.0.\n",
      "  dantzig solved it with objective 125047848.0.\n",
      "  gurobi solved it with objective 125047848.0.\n",
      "n=30\n",
      "  add_circuit solved it with objective 95776200.0.\n",
      "  miller_tucker_zemlin solved it with objective 95776200.0.\n",
      "  dantzig solved it with objective 95776200.0.\n",
      "  gurobi solved it with objective 95776200.00000001.\n",
      "n=40\n",
      "  add_circuit solved it with objective 90770264.0.\n",
      "  miller_tucker_zemlin solved it with objective 90770264.0.\n",
      "  dantzig solved it with objective 90770264.0.\n",
      "  gurobi solved it with objective 90770264.0.\n",
      "n=50\n",
      "  add_circuit solved it with objective 87023666.0.\n",
      "  gurobi solved it with objective 87023666.0.\n",
      "  miller_tucker_zemlin dropped out.\n",
      "  dantzig dropped out.\n",
      "n=60\n",
      "  add_circuit solved it with objective 86011310.0.\n",
      "  gurobi solved it with objective 86011310.0.\n",
      "n=70\n",
      "  add_circuit solved it with objective 83601078.0.\n",
      "  gurobi solved it with objective 83601078.0.\n",
      "n=80\n",
      "  gurobi solved it with objective 88234698.0.\n",
      "  add_circuit dropped out.\n",
      "n=90\n",
      "  gurobi solved it with objective 80044620.0.\n",
      "n=100\n",
      "  gurobi solved it with objective 70987602.0.\n",
      "n=110\n",
      "  gurobi solved it with objective 72344612.0.\n",
      "n=120\n",
      "  gurobi solved it with objective 71059516.0.\n",
      "n=130\n",
      "  gurobi solved it with objective 76303294.0.\n",
      "n=140\n",
      "  gurobi solved it with objective 78506514.0.\n",
      "n=150\n",
      "  gurobi solved it with objective 78223603.90126976.\n",
      "n=160\n",
      "  gurobi solved it with objective 72123102.00000006.\n",
      "n=170\n",
      "  gurobi solved it with objective 72509394.0.\n",
      "n=180\n",
      "  gurobi dropped out.\n"
     ]
    }
   ],
   "source": [
    "challengers = {\n",
    "    \"add_circuit\": solve_with_add_circuit,\n",
    "    \"miller_tucker_zemlin\": solve_with_miller_tucker_zemlin,\n",
    "    \"dantzig\": solve_with_dantzig,\n",
    "    \"gurobi\": solve_with_gurobi,\n",
    "}\n",
    "n = STEP_SIZE\n",
    "while challengers:\n",
    "    print(f\"n={n}\")\n",
    "    G, points = generate_random_euclidean_graph(n)\n",
    "    dropouts = []\n",
    "    for name, solver in challengers.items():\n",
    "        solved, obj = solver(G)\n",
    "        if not solved:\n",
    "            dropouts.append(name)\n",
    "        else:\n",
    "            print(f\"  {name} solved it with objective {obj}.\")\n",
    "    for name in dropouts:\n",
    "        print(f\"  {name} dropped out.\")\n",
    "        del challengers[name]\n",
    "    n += STEP_SIZE"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Benchmark\n",
    "\n",
    "Now let us do a slightly more serious benchmark, by running each model multiple times and averaging the results.\n",
    "This is still not a very good benchmark, because the dropout can skew the results.\n",
    "It would be better to let each approach try to solve every instance, and then compare the results by looking at how many instances of each size could still be solved.\n",
    "However, our approach here is more time efficient."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "{'dantzig': 48.0, 'miller_tucker_zemlin': 48.0, 'add_circuit': 110.0, 'gurobi': 225.0}\n"
     ]
    }
   ],
   "source": [
    "def run_dropout_benchmark():\n",
    "    challengers = {\n",
    "        \"add_circuit\": solve_with_add_circuit,\n",
    "        \"miller_tucker_zemlin\": solve_with_miller_tucker_zemlin,\n",
    "        \"dantzig\": solve_with_dantzig,\n",
    "        \"gurobi\": solve_with_gurobi,\n",
    "    }\n",
    "    dropout_points = dict()\n",
    "    n = STEP_SIZE\n",
    "    while challengers:\n",
    "        G, points = generate_random_euclidean_graph(n)\n",
    "        dropouts = []\n",
    "        for name, solver in challengers.items():\n",
    "            solved, obj = solver(G)\n",
    "            if not solved:\n",
    "                dropouts.append(name)\n",
    "        for name in dropouts:\n",
    "            del challengers[name]\n",
    "            dropout_points[name] = n\n",
    "        n += STEP_SIZE\n",
    "    return dropout_points\n",
    "\n",
    "\n",
    "dropout_points = [run_dropout_benchmark() for _ in range(10)]\n",
    "# compute the average\n",
    "dropout_points = {\n",
    "    name: sum(d[name] for d in dropout_points) / len(dropout_points)\n",
    "    for name in dropout_points[0]\n",
    "}\n",
    "print(dropout_points)"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "mo310",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.12.3"
  },
  "orig_nbformat": 4
 },
 "nbformat": 4,
 "nbformat_minor": 2
}
